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Abstract 

We show that the self-gravitating gas at thermal equilibrium has an infinite 
volume limit in the three ensembles (GCE, CE, MCE) when (A^, V) oo, 
keeping N/V^^^ fixed, that is, with r] = ^f/kY fixed. We develop Monte Carlo 
simulations, analytic mean field methods (MF) and low density expansions. 
We compute the equation of state and find it to be locally p{r) = T pv{r), 
that is a local ideal gas equation of state. The system is in a gaseous phase 
for Tj < TjT = 1.51024. . . and collapses into a very dense object for rj > rjT 
in the CE with the pressure becoming large and negative. The isothermal 
compressibility diverges ed, r] = tjt- We compute the fiuctuations around 
mean field for the three ensembles. We show that the particle distribution 
can be described by a Haussdorf dimension 1 < D < 3. 

Nous montrons que la limite de volume infini existe pour le gaz auto- 
gravitant a I'equilibre thermique dans les trois ensembles (EGC,EC,EMC) 
quand (iV, V) — >• oo, avec N/V^^^ fixe, c'est a dire r] = ^t/sy Nous util- 
isons les simulations Monte Carlo, la methode du champ moyen et les devel- 
oppements a basse densite. Nous calculous I'equation d'etat et nous trouvons 
qu'elle est localement p{r) = T pv{r), c'est dire, I'quation d'un gaz parfait 
local. Le system est dans une phase gazeuse pour f] < rjT = 1.51024 . . . 
et s'effondre dans un objet tres dense pour f] > r]T dans I'ensemble canon- 
ique avec une pression grande et negative. La compressibilite isothermique 
diverge k rj = rjT- Nous calculous les fluctuations autour du champ moyen 
pour les trois ensembles. Nous montrons que la distribution des particules 
est decrite par une dimension de Haussdorf \ < D < 3. 
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Physical systems at thermal equilibrium are usually homogeneous. This 
is the case for gases with short range intermolecular forces (and in absence 
of external fields). In such cases the entropy is maximum when the system 
homogenizes. 

When long range interactions as the gravitational force are present, even 
the ground state is inhomogeneous. In this case, each element of the sub- 
stance is acted on by very strong forces due to distant particles in the gas. 
Hence, regions near to and far from the boundary of the volume occupied by 
the gas will be in very different conditions, and, as a result, the homogeneity 
of the gas is destroyed . The state of maximal entropy for gravitational sys- 
tems is inhomogeneous. This basic inhomogeneity suggested us that fractal 
structures can arise in a self-interacting gravitational gas [21 El IH El El • 

The inhomogeneous character of the ground state for gravitational sys- 
tems explains why the universe is not going towards a 'thermal death'. A 
'thermal death' would mean that the universe evolves towards more and more 
homogeneity. This can only happen if the entropy is maximal for an homoge- 
neous state. Instead, it is the opposite what happens, structures are formed 
in the universe through the action of the gravitational forces as time evolves. 

Usual theorems in statistical mechanics break down for inhomogeneous 
ground states. For example, the specific heat may be negative in the mi- 
crocanonical ensemble (not in the canonical ensemble where it is always 
positive) |T]. 

As is known, the thermodynamic limit for self-gravitating systems does 
not exist in its usual form (A^ ^ oo, V ^ oo, N/V = fixed). The system 
collapses into a very dense phase which is determined by the short distance 
(non-gravitational) forces between the particles. However, the thermody- 
namic functions exist in the dilute limit |Sl IH El 

A^ 

N oo , V oo , = fixed , 

where V stands for the volume of the box containing the gas. In such a limit, 
the energy E, the free energy and the entropy turns to be extensive. That 
is, we find that they take the form of A^ times a function of the intensive 
dimensionless variables: 

Gm^N ^ EL 

V = — — '^^ ^ 



LT ^ Gm^N'^ 

where rj and ^ are intensive variables. Namely, rj and ^ stay finite when A^ 
and V = tend to infinite. The variable rj is appropriate for the canonical 
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ensemble and ^ for the microcanonical ensemble. Physical magnitudes as 
the specific heat, speed of sound, chemical potential and compressibility only 
depend on 1] or ^. The variables r] and ^, as well as the ratio N/L, are 
therefore intensive magnitudes. The energy, the free energy, the Gibbs free 
energy and the entropy are of the form N times a function of rj. These 
functions of rj have a finite N = oo limit for fixed 1] (once the ideal gas 
contributions are subtracted). Moreover, the dependence on rj in all these 
magnitudes express through a single universal function f(^)- The variable r) 
is the ratio of the characteristic gravitational energy ^ and the kinetic 
energy T of a particle in the gas. For 77 = the ideal gas is recovered. 

In refs. j3| HI and [Uj we have thoroughly studied the statistical mechan- 
ics of the self-gravitating gas. That is, our starting point is the partition 
function for non-relativistic particles interacting through their gravitational 
attraction in thermal equilibrium. We study the self-gravitating gas in the 
three ensembles: microcanonical (MCE), canonical (CE) and grand canonical 
(GCE). We performed calculations by three methods: 

• By expanding the partition function through direct calculation in pow- 
ers of 1/^ and 1] for the MCE and CE, respectively. These expressions 
apply in the dilute regime (^ ^ 1 , <^ 1) and become identical for 
both ensembles for 00. At rj = = 1/^ we recover the ideal gas 
behaviour. 

• By performing Monte Carlo simulations both in the MCE and in the 
CE. We found in this way that the self-gravitating gas collapses at 
a critical point which depends on the ensemble considered. As shown 
in fig. Qthe collapse occurs first in the canonical ensemble (point T). 
The microcanonical ensemble exhibits a larger region of stability that 
ends at the point MC (fig. H}. Notice that the physical magnitudes 
are identical in the common region of validity of both ensembles within 
the statistical error. Beyond the critical point T the system becomes 
suddenly extremely compact with a large negative pressure in the CE. 
Beyond the point MC in the MCE the pressure and the temperature 
increase suddenly and the gas collapses. The phase transitions at 
T and at MC are of zeroth order since the Gibbs free energy has 
discontinuities in both cases. 

• By using the mean field approach we evaluate the partition function 
for large A^. We do this computation in the grand canonical, canonical 
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Figure 1: f{f]) = PV/[NT] as a function of rj by Monte Carlo simulations 
for the microcanonical and canonical ensembles (A^ = 2000). Both curves 
coincide within the statistical error till the point T. 
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and microcanonical ensembles. In the three cases, the partition func- 
tion is expressed as a functional integral over a statistical weight which 
depends on the (continuous) particle density. These statistical weights 
are of the form of the exponential of an 'effective action' proportional 
to A^. Therefore, the N —>■ oo limit follows by the saddle point method. 
The saddle point is a space dependent mean field showing the inhomo- 
geneous character of the ground state. Corrections to the mean field 
are of the order and can be safely ignored for N ^ 1 except near 
the critical points. These mean field results turned out to be in excel- 
lent agreement with the Monte Carlo results and with the low density 
expansion. 

We calculate the saddle point (mean field) for spherical symmetry and we ob- 
tain from it the various physical magnitudes (pressure, energy, entropy, free 
energy, specific heats, compressibilities, speed of sound and particle density). 
Furthermore, we computed in ref. |3] the determinants of small fluctua- 
tions around the saddle point solution for spherical symmetry for the three 
statistical ensembles. 

When any small fluctuation around the saddle point decreases the statis- 
tical weight in the functional integral, the saddle point dominates the func- 
tional integral and the mean field approach can be valid. In that case, the 
determinant of small fluctuations is positive. A negative determinant of small 
fluctuations indicates that some fluctuations around the saddle point are in- 
creasing the statistical weight in the functional integral and hence the saddle 
point does not dominate the partition function. The mean field approach 
cannot be used when the determinant of small fluctuations is negative. We 
find analytically in the CE that the determinant of small fluctuations van- 
ishes at the point rjc = 1.561764 and becames negative a.t f] > rjc 0- (The 
point 7]c is indicated C in Figure 1). We find that the CE specific heat cy 
at the point C diverges as j^j: 

cv ±0.63572 ...{vc- vY^'^ - 0.19924 . . . + 0(Vr/c - rj) 

The (+) sign refers to the positive (first) branch, and the (-) sign to the 
negative (second) branch (between the points C and MC). 

However, it must be noticed that the instability point is located at 
T) = r)T < rjc, as shown by both, mean field and Monte Carlo computations. 
(The point rjT is indicated T in Figure 1). The onset of instability in the 
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canonical ensemble coincides with the point where the isothermal compress- 
ibility diverges. The isothermal compressibility is positive from rj = 
till rj = TjT = 1.51024 .... At this point as well as the specific heat at 
constant pressure cp diverge and change their signs Moreover, at this 
point the speed of sound at the center of the sphere becomes imaginary |1] . 
Therefore, small density fiuctuations will grow exponentially in time instead 
of exhibiting oscillatory propagation. Such a behaviour leads to the collapse 
of the gas into a extremely compact object. Monte-Carlo simulations confirm 
the presence of this instability at rjT = 1.510 ... in the canonical ensemble 
and the formation of the collapsed object 

The collapse in the Grand Canonical ensemble (GCE) occurs for a smaller 
value of 7] = rjGc = 0.49465 . . ., while in the micro canonical ensemble, the 
collapse arrives later, in the second sheet, at 77 = rjMC = 1.25984. . .. 

We find that the Monte Carlo simulations for self-gravitating gas in the 
CE and the MCE confirm the stability results obtained from mean field. 

The saddle point solution is identical for the three statistical ensembles. 
This is not the case for the fluctuations around it. The presence of constraints 
in the CE (on the number of particles) and in the MCE (on the energy and 
the number of particles) changes the functional integral over the quadratic 
fluctuations with respect to the GCE. 

The saddle point of the partition function turns out to coincide with the 
hydrostatic treatment of the self-gravitating gas [7j (which is usually known 
as the 'isothermal sphere' in the spherically symmetric case). 

We flnd that the Monte Carlo simulations (describing thermal equi- 
librium) are much more efficient than the A^-body simulations integrating 
Newton's equations of motion. (Indeed, the integration of Newton's equa- 
tions provides much more detailed information than the one needed in ther- 
mal equilibrium investigations). Actually, a few hundreds of particles are 
enough to get quite accurate results in the Monte Carlo simulations (except 
near the collapse points). Moreover, the Monte Carlo results turns to be 
in excellent agreement with the mean fleld calculations up to very small 
corrections of the order (1/A^). Our Monte Carlo simulations are performed 
in a cubic geometry. The equilibrium conflgurations obtained in this manner 
can thus be called the 'isothermal cube'. 

In summary, the picture we get from our calculations using these three 
methods show that the self-gravitating gas behaves as a perfect gas for 
?7 — 0, 1/.^ ^ 0. When rj and 1/^ grow, the gas becomes denser till it 
suddenly condenses into a high density object at a critical point GC, C or 
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Microcanonical, collapsed, zi=-0.6 




Figure 2: Average particle distribution in the collapsed phase from Monte 
Carlo simulations with 2000 particles in the microcanonical ensemble for 
^ = -0.6, T] = 0.43, pV/[NT] = 0.414. 

MC depending upon the statistical ensemble chosen. 

T] is related with the Jeans' length dj of the gas through r] = 3 {L/dj^. 
Hence, when t] goes beyond t]t, the length of the system becomes larger than 
dj / \J TfT /3. The collapse at T in the CE is therefore a manifestation of the 
Jeans' instability. 

In the MCE, the determinant of fluctuations vanishes at the point MC. 
The physical states beyond MC are collapsed configurations as shown by the 
Monte Carlo simulations (see fig. El). Actually, the gas collapses in the Monte 
Carlo simulations slightly before the mean field prediction for the point MC. 
The phase transition at the microcanonical critical point MC is the so called 
gravothermal catastrophe |H]. 

The gravitational interaction being attractive without lower bound, a 
short distance cut-off {A) must be introduced in order to give a meaning to 
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the partition function. We take the gravitational force between particles as 
—G IV? /r"^ for r > A and zero for r < A where r is the distance between the 
two particles. We show that the cut-off effects are negligible in the N = oo 
limit. That is, once we set = cxd with fixed rj, all physical quantities are 
finite in the zero cut-off limit {A = 0). The cut-off effects are of the order 
A^ / L"^ and can be safely ignored. 

In ref. j2] we expressed all global physical quantities in terms of a sin- 
gle function /(r/). Besides computing numerically f{ri) in the mean field 
approach, we showed that this function obeys a first order non-linear differ- 
ential equation of first Abel's type. We obtained analytic results about firf) 
from the Abel's equation, /(r^) exhibits a square-root cut at rjc-, the 
critical point in the CE. The first Riemann sheet is realized both in the CE 
and the MCE, whereas the second Riemann sheet (where cy < 0) is only 
realized in the MCE. f{7]) has infinitely many branches in the t] plane but 
only the first two branches are physically realized. Beyond MC the states 
described by the mean field saddle point are unstable. 

We plot and analyze the equation of state, the energy, the entropy, the 
free energy, cy and the isothermal compressibility (figs. 9-13 in ref. [HI)- Most 
of these physical magnitudes were not previously computed in the literature 
as functions of t]. 

We find analytically the behaviour of f{ri) near the point rjc in mean 
field, 

fiv) ''=^ I + 0.27137 . . . + 0.27763 . . . ir]c - v) + O [{vc - vf^' ■ 

This exhibits a square root branch point singularity at r]c- This shows that 
the specific heat at constant volume diverges at rjc as {rjc — V ^ Vc- 

However, it must be noticed, that the specific heat at constant pressure and 
the isothermal compressibility? both diverge at the point rjT as {rjT — ri)~^ . 
These mean field results apply for \ri — ric\ <^ 1 <^ N\ri — ric\- Fluctuations 
around mean field can be neglected in such a regime. 

The Monte Carlo calculations permit us to obtain /(?]) in the collapsed 
phase. Such result cannot be obtained in the mean field approach. The mean 
field approach only provides information as /(//) in the dilute gas phase. 

For the self-gravitating gas, we find that the Gibbs free energy $ is not 
equal to times the chemical potential and that the thermodynamic poten- 
tial Q is not equal to —PV as usual ^J. This is a consequence of the dilute 
thermodynamic limit A^ oo, L oo, N/L =fixed. 
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We computed the determinant of small fluctuations around the saddle 
point solution for spherical symmetry in all three statistical ensembles. In 
the spherically symmetric case, the determinant of small fluctuations is writ- 
ten as an infinite product over partial waves. The S and P wave determinants 
are written in closed form in terms of the saddle solution. The determinants 
for higher partial waves are computed numerically. All partial wave deter- 
minants are positive definite except for the S-wave P]. The reason why 
the fluctuations are different in the three ensembles is rather simple. The 
more contraints are imposed the smaller becomes the space of fluctuations. 
Therefore, in the grand canonical ensemble (GCE) the system is more free 
to fluctuate and the phase transition takes place earlier than in the micro- 
canonical (MCE) and canonical ensembles (CE). For the same reason, the 
transition takes place earlier in the CE than in the MCE. 

The conclusion being that the mean field correctly gives an excellent 
description of the thermodynamic limit except near the critical points (where 
the small fluctuations determinant vanishes); the mean field is valid for N\ri — 
Vcritl ^ 1- The vicinity of the critical point should be studied in a double 
scaling limit oo, 77 —> rjcru- Critical exponents are reported in ref. 

[2| for t] ^ rjc using the mean field. These mean field results apply for 
\v ~ Vc\ 1 A^l'? — ''7c| with N ^ 1. Fluctuations around mean field can 
be neglected in such a regime. 

We computed local properties of the gas in That is, the local energy 
density e(r), local particle density, local pressure and the local speed of sound. 
Furthermore, we analyze the scaling behaviour of the particle distribution 
and its fractal (Haussdorf) dimension i4j. 

The particle distribution pvio) proves to be inhomogeneous (except for 
?7 <^ 1) and described by an universal function of rj, the geometry and the 
ratio f = q/R, R being the radial size. Both Monte Carlo simulations and 
the Mean Field approach show that the system is inhomogeneous forming a 
clump of size smaller than the box of volume V [see fig. El here and figs. 3, 
5 and 6 in [3]]. 

The particle density in the bulk behaves as pv{(f) — r^~^. That is, the 
mass M{R) enclosed on a region of size R vary approximately as 

M{R) R^ . 

D slowly decreases from the value D = 3 for the ideal gas (77 = 0) till 
D = 0.98 in the extreme limit of the MC point, D takes the value 1.6 at 
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rjc, [see Table 1]. This indicates the presence of a fractal distribution with 
Haussdorf dimension D. 

Our study of the statistical mechanics of a self-gravitating system indi- 
cates that gravity provides a dynamical mechanism to produce fractal struc- 
tures [mill. 

The average distance between particles monotonically decrease with rj in 
the first sheet. The mean field and Monte Carlo are very close in the gaseous 
phase whereas the Monte Carlo simulations exhibit a spectacular drop in the 
average particle distance at the clumping transition point T. In the second 
sheet (only described by the MCE) the average particle distance increases 
with T] 

We find that the local equation of state is given by 

p(f) =Tpvir} . 

We have derived the equation of state for the self-gravitating gas. It is 
locally the ideal gas equation, but the self-gravitating gas being inhomo- 
geneous, the pressure at the surface of a given volume is not equal to the 
temperature times the average density of particles in the volume. In partic- 
ular, for the whole volume: PV/[NT] = f{r]) < 1 (the equality holds only 
for ?7 = 0). 

Notice that we have found the local ideal gas equation of state p{r) = 
T pv{r) for purely gravitational interaction between particles. Therefore, 
equations of state different from this one, (as often assumed and used in 
the literature for the self-gravitating gas), necessarily imply the presence of 
additional non-gravitational forces. 

The local energy density e(r) turns out to be an increasing function of r 
in the spherically symmetric case. The energy density is always positive on 
the surface, whereas it is positive at the center for < < r/s = 1.07783 . . ., 
and negative beyond the point rj = rj^ = 1.07783 .... 

The local speed of sound v'^{r) is computed in the mean field approach 
as a function of the position for spherical symmetry and long wavelengths. 
Vg{r) diverges at t] = rjT = 1.51024. . . in the first Riemann sheet. Just be- 
yond this point Vg{r) is large and negative in the bulk showing the strongly 
unstable behaviour of the gas for such range of values of rj. 

Moreover, we have shown the equivalence between the statistical me- 
chanical treatment in the mean field approach and the hydrostatic description 
of the self-gravitating gas [7]. 
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The success of the hydrodynamical description depends on the value of 
the mean free path (/) compared with the relevant sizes in the system. / must 
be -C 1. We compute the ratio l/a (Knudsen number), where a is a length 
scale that stays fixed for N ^ oo and show that Ija ^ . This result 
ensures the accuracy of the hydrodynamical description for large N . 

Furthermore, we have computed in ref. j3] several physical magnitudes 
as functions of r] and r which were not previously computed in the literature 
as the speed of sound, the energy density, the average distance between 
particles and we notice the presence of a Haussdorf dimension in the particle 
distribution. 

The statistical mechanics of a selfgravitating gas formed by particles with 
different masses is thoroughly investigated in ref. [H] while the selfgravitating 
gas in the presence of the cosmological constant is thoroughly investigated 
in refs. [H]. In ref.^U] the Mayer expansion for the selfgravitating gas is 
investigated in connection with the stability of the gaseous phase. 
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TABLE 1. The Fractal Dimension D and the proportionality coefficient 
C as a function of rj from a fit to the mean field results according to M(r) ~ 
C r^. The point rjMc is on the second Riemann sheet. 
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